Instantaneous Gelation in Smoluchowski's Coagulation Equation Revisited 
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We study the solutions of the Smoluchowski coagulation equation with a regularisation term 
which removes clusters from the system when their mass exceeds a specified cut-off size, M. We 
focus primarily on collision kernels which would exhibit an instantaneous gelation transition in 
the absence of any regularisation. Numerical simulations demonstrate that for such kernels with 
monodisperse initial data, the regularised gelation time decreases as M increases, consistent with 
the expectation that the gelation time is zero in the unregularised system. This decrease appears to 
be a logarithmically slow function of M , indicating that instantaneously gelling kernels may still be 
justifiable as physical models despite the fact that they are highly singular in the absence of a cut-off. 
We also study the case when a source of monomers is introduced in the regularised system. In this 
case a stationary state is reached. We present a complete analytic description of this regularised 
stationary state for the model kernel, K(m\,m2) = max {mi, 7712}", which gels instantaneously 
when M — >■ 00 if > 1. The stationary cluster size distribution decays as a stretched exponential 
for small cluster sizes and crosses over to a power law decay with exponent i/ for large cluster sizes. 
The total particle density in the stationary state slowly vanishes as (log(Af''~^))~^''^ when M — >■ 00. 
The approach to the stationary state is non-trivial : oscillations about the stationary state emerge 
from the interplay between the monomer injection and the cut-off, M , which decay very slowly when 
M is large. A quantitative analysis of these oscillations is provided for the addition model which 
describes the situation in which clusters can only grow by absorbing monomers. 

PACS numbers: 83.80.Jx 



I. INTRODUCTION 

The determination of the statistical time evolution of 
an ensemble of particles undergoing irreversible binary 
coagulation is a problem which arises in many different 
contexts in the physical, chemical and biological sciences. 
See ^ for a list of applications. Of pre-eminent interest 
is the cluster size distribution, Nm (t) , which specifies the 
average density of particles of mass, m, at a given time 
t. Often one aims to derive Nm (t) from a given model of 
the microscopic dynamics governing cluster coagulation. 
This problem has been extensively studied for almost a 
century starting with the seminal work of Smoluchowski 
[2] who showed that, at the mean field level, the clus- 
ter size distribution of a statistically homogeneous sys- 
tem evolves according to the Smoluchowski coagulation 
equation: 

dtNm^ / dmiK{mi,m- mi)NmiNm-m^ (1) 

f°° J 

— 2Nm I dmiK{m,mi)N,n H 6{m~ itiq). 

Jq mo 

Here the coagulation kernel, K (mi, 1112), is proportional 
to the probability rate for a cluster of mass mi merging 
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with a cluster of mass 7712. At this level of description, 
it encodes all relevant information about the underlying 
micro-physics. J is the rate of injection of monomers, 
having mass toq, into the system. J may be zero de- 
pending on the application. In this paper, we take Eq. ([T]) 
as our departure point. It is expected to apply when the 
clusters remain well-mixed (despite aggregation) , but the 
elucidation of the exact conditions under which Eq. ([T]) is 
rigorously obtained as the mean field limit of an under- 
lying stochastic process can be a subtle question. Math- 
ematically inclined readers are referred to the review by 
Aldous [3] for some common mathematical perspectives 
on this issue. A physical example of a situation in which 
the applicability of Eq. ([!]) fails due to the generation of 
spatial correlations between particles by diffusive fluctu- 
ations is discussed in detail in [3]. 

In many applications the microphysics is scale invari- 
ant, at least over some range of cluster sizes. That is 
to say, the coagulation kernel is a homogeneous function 
of its arguments, the degree of homogeneity of which we 
shall denote by A: 

K{cmi,cm2) — c^K{mi,ra2). (2) 

In many cases, such kernels result in solutions of Eq. ([T]) 
which exhibit self-similarity. This means that the cluster 
size distribution asymptotically takes the form N{m, t) ~ 
s{t)°'F{m/ s{t)) where s{t) is a characteristic cluster size 
which grows in time, a is a dynamical scaling exponent 
and ~ denotes the scaling limit: t — > 00 and to — > cxi 
with m/s(t) fixed. Furthermore, it has been proven for 
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FIG. 1. (Color online) Total mass density, A4i{t), for 
K{mi,m2) = (77117712)"'"''"'. 



FIG. 2. (Color online) Total mass density, A4i{t), for 



the constant, sum and product kernels that such self- 
similar solutions are attracting for any initial data "with 
finite (A-l- 1)*'' moment [Sj. Much "work in the physics lit- 
erature on the theory of coagulation, follo"wing the work 
of Van Dongen and Ernst , has focused on the prob- 
lem of determining the properties of the scaling function 
F(x) and the exponent a from the scaling properties of 
K {mi, 1712)- An almost complete scaling theory of Eq. ([T]) 
is now known. See the review by Leyvraz [7] for a modern 
discussion. 

A key feature of this scaling theory, and one which is of 
considerable importance for what follows, is the fact that 
the scaling properties of Eq. ([T]) are often very sensitive 
to the rate of coagulation of clusters of widely different 
masses. In order to parameterise this dependence in a 
general way, and following the notation of [7], we intro- 
duce the scaling exponents fi and v which specify the 
asymptotic behaviour of the coagulation kernel in the 
limit where the mass of one cluster greatly exceeds that 
of the other: 

K {mi, 1112) ^ m'^ 7712 'tni<tJ7i2- (3) 

Clearly from Eq. ([2]), we must have 11 + ly ~ \. We shall 
use the notation A4a{t) to denote the a-moment of the 
cluster size distribution: 

Ma{t)^ / m°'N„r{t)dm. (4) 

The first moment, the total mass density, A4i{t) — 
/q m Nm{t) dm, is governed by conservation of material 
and is thus formally conserved by Eq. (IT]) when J — 0. 
If J 7^ we simply have Mi{t) = Mi{oY+ Jt. This ac- 
cords with the intuition imparted by the mass-conserving 
character of the individual coagulation events. This intu- 
ition is challenged, however, by one of the more interest- 
ing phenomena to emerge from the scaling theory of the 



Smoluchowski equation; the so-called gelation transition. 
Originally conjectured by Lushnikov [5] and Ziff [|] and 
put on a firmer theoretical footing by Van Dongen and 
Ernst [ilO], gelation refers to the fact that, for kernels 
having A > 1, mass-conservation can be spontaneously 
broken in finite time. That is to say, there exists a time 
t* , known as the gelation time such that 

Mi{t) < I mN^{0)dm t > t* . 
Jo 

The characteristic cluster size, s{t), typically diverges at 
the gelation time. The "missing" mass can be interpreted 
as going into a cluster of infinite mass or "gel" . The gen- 
eration of arbitrarily large clusters in finite time might 
seem surprising but can be perfectly physical in some 
cases. One example is polymer gelation in which clus- 
ters merge by the formation of crosslinks and therefore 
do not need to move in order to coalesce. For other ex- 
amples of kernels having A > 1, the solution of Eq. ([T]) 
describes only the intermediate asymptotics of the un- 
derlying physical system over some range of cluster sizes 
and for times less than the gelation time. Once this inter- 
mediate asymptotic range has been identified, there is no 
conceptual problem with the loss of mass conservation in 
the Smoluchowski or the generation of infinite clusters. 
Indeed, the gelation transition has even been observed 
experimentally, for example in polymer aggregation JJJ, 
and found to exhibit dynamics in reasonable agreement 
with the predictions of the Smoluchowski equation. 

It turns out, however, that this is not the full story. 
Soon after the discovery of the gelation phenomenon, it 
was realised by Hendriks, Ernst and Ziff [T^ that the ker- 
nel K{mi,m2) = (miTO2)^ produced a series expansion 
for the second moment of the cluster size distribution 
which seemed to have zero radius of convergence in time. 
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This led the authors to suggest the possibilty of gelation 
occuring at time 0+. Detailed study of the scaling prop- 
erties of the Smoluchowski equation subsequently led to 
formal arguments |13], later made rigourous |14j . which 
show that for coagulation kernels having exponent v > 1, 
the gelation time is actually zero. This surprising phe- 
nomenon is referred to as instantaneous gelation. Even 
more surprisingly, the gelation process can be complete in 
the sense that Aii{t) = for i > 0. See [TS] for a mathe- 
matical discussion of this point. A situation in which all 
mass vanishes from the system in time 0-1- clearly can- 
not describe even the intermediate asymptotics of any 
physical coagulation problem. For applications such as 
the coagulation of polymers or colloidal aggregates, the 
fact that the available surface area of an aggregate can- 
not grow faster than its volume means that the exponent 
J/ cannot be greater than 1. If the only applications of 
the Smoluchowski equation came from polymer science, 
the phenomenon of instantaneous gelation would be re- 
garded as a mathematical pathology which need not con- 
cern physicists. Yet there are models for which it can 
plausibly be argued that the exponent v is greater than 
1. One such example is the astrophysical phenomenon of 
gravitational clustering HHITT] which is thought to play a 
role in determining the large scale matter distribution of 
the universe. A second important example is that of dif- 
ferential sedimentation of water droplets falling at their 
terminal velocity [IHl [H], one of the processes respon- 
sible for the observed droplet size distribution in clouds 
[201 121j. Furthermore, there are even proposed heuris- 
tic solutions of Eq. ([!]) in the literature for such models 
[17L 119] which seem reasonably supported by numerical 
simulations. The question of how this is possible, given 
the known mathematical results on instantaneous gela- 
tion discussed above, is the principal topic of this paper. 
It should be clear from the outset that the presence of an 
instantaneous gelation transition for a particular kernel 
indicates that the underlying physical model omits some 
process which becomes important for short timescales or 
for large masses. For example, in the case of droplet 
coalescence in clouds, the Stokes assumption for the ter- 
minal velocity of a droplet is ultimately responsible for 
the exponent v exceeding 1. This assumption ceases to 
be valid for sufficiently large droplets (see [18] and the 
references therein). 



II. INSTANTANEOUS GELATION IN THE 
REGULARISED SYSTEM 

The gelation phenomenon in particle systems is best 
understood by considering a regularisation of the system 
and studying the behaviour as this regularisation is re- 
moved. Two natural regularisations can be found in the 
literature. One approach is to consider the stochastic dy- 
namics of a finite number of particles as has been done by 
Lushnikov [52] for the product kernel. Another approach 
is to introduce a mass cut-off. A/, into the Smoluchowski 



equation. This has been done by Filbet and Laurengot 
[23j . These are different regularisations but both show 
singular behaviour as the regularisation is removed when 
A > 1. 

In this paper, we regularise by the latter method. As in 
[53] , the cut-off is introduced in such a way that clusters 
having mass m > M are removed from the system: 

1 

dtN„i = - J dmiK{mi,m~ mi)N^^N^_jn, (5) 

/M—m 
d'miK{m, mi)Nmi + J S{m — 1) 

- Dm [N„,{t)] 

where 

Dm [N,n {t)]^N,n dm,K{m,m,)N^^ (6) 

J M-m 

describes the removal of clusters having m > M. This 
regularisation explicitly breaks mass conservation. Fur- 
thermore, we shall implicitly measure all masses in terms 
of the monomer mass from now on so that the lower cut- 
off, mg, is set equal to 1. Clearly this is not the only 
way in which the problem can be regularised. In particu- 
lar, by omitting Dm [Nm{t)] from Eq. ^ one would ob- 
tain an explicitly conservative regularisation. One would 
expect the behaviour to be completely different in this 
case. The differences between conservative and non- 
conservative regularisations for a related aggregation- 
fragmentation model arising in the kinetics of waves are 
explored in [53]. The fact that the gelation transition 
involves a loss of mass from the system suggests that the 
non-conservative regularisation is more natural. While 
we implement this regularisation as a hard cut-off, M, 
we would expect the behaviour described below to be 
qualitatively the same for smoother regularisations pro- 
vided that the non-conservative property is retained. 

Figs. [T] and [2| show the mass contained in this 
regularised system as a function of time for a se- 
quence of values of M for the non-gelling kernel 
K{mi,m2) = (mim2)^/^ (Fig. Ill) and for the gelling 
kernel K{mi,m,2) = (mim2)'^/'* (Fig. J2|. These were 
obtained by numerical solution of Eq. ([sT with monodis- 
perse initial data. 

All numerical results in this paper were obtained us- 
ing the algorithm developed by Lee in [53 [5S]. This 
algorithm involves coarse-graining the masses into bins 
whose width increases exponentially with mass and then 
using the Smoluchowski equation to approximate the 
mass transferred per timestep as a result of the aggre- 
gation of clusters in each pair of bins. The approxima- 
tion explicitly enforces mass conservation. Full details of 
the coarse-graining and formulae for the computation of 
mass transfer are provided in |25| . In a slight modifica- 
tion of Lee's algorithm, the regularisation described by 
Eq. (p| is implemented by introducing an infinite width 
bin, pVf , oo] , which accumulates any mass transferred to 
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FIG. 3. (Color online) Total mass density, A4i{t), for the 
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clusters larger than M. We also modified the timestep- 
ping procedure, replacing Lee's original explicit timestep- 
ping method with an adaptive implicit trapezoidal rule. 
This is essential because the Smoluchowski equation be- 
comes increasingly stiff for larger values of A to the ex- 
tent that the results presented here could not be obtained 
with explicit timestepping. The implicit trapezoidal rule 
requires the solution of a set of nonlinear equations at 
each timestep which was done using the GSL im- 
plementation of the Rosenbrock algorithm [5H]. Explicit 
formulae and further discussion of these modifications of 
Lee's algorithm can be found in [24] . 

We see that for the non-gelling system, mass conser- 
vation is restored in the limit M ^ oo whereas for the 
gelling system, it is not. This latter situation will be 
recognisable to readers familiar with the theory of turbu- 
lence where it is widely believed that energy conservation 
is not restored when the limit of zero viscosity is taken in 
the Navier-Stokes equations [29) . a phenomenon referred 
to as the dissipative anomaly. A similar phenomenon is 
observed in the kinetics of wave turbulence ^3]. There 
is no physical contradiction in the fact that mass con- 
servation is broken in the Smoluchowski equation in the 
gelling regime. It simply means that for times larger than 
the gelation time, the underlying conservative coagula- 
tion dynamics must be modified for the largest clusters. 

Let us now consider what happens when an instan- 
taneously gelling kernel is inserted into the regularised 
Smoluchowski equation. An archetypal instantaneously 
gelling kernel, which we shall use extensively in what fol- 
lows, is the generalised sum kernel: 

i^,(mi,TO2) = m}+' -Km^+^ (7) 

For this kernel \ = v ^ 1 + e and /i = 0. According to the 
classification of Van Dongen and Ernst, it is nongelling 



for e < and instantaneously gelling for e > 0. In the 
marginal case, e = 0, it is the simple sum kernel which 
is exactly solvable, at least in the absence of a source 
of mononers [7 , and turns out to be non-gelling. Fig. [s] 
shows the sol mass in the regularised system with the gen- 
eralised sum kernel with — \ for a sequence of increas- 
ing values of the regularisation mass, M . In the presence 
of the cut-off, the regularized gelation time, t\,j > 0, is 
clearly identifiable. This regularised gelation time, al- 
though finite, decreases as M^ax increases. Extrapolat- 
ing the behaviour seen in Fig. [3] it is plausible that the 
gelation time vanishes as M — oo consistent with the 
expectation that the unregularised system exhibits com- 
plete instantaneous gelation. Instantaneous gelation has 
not, to the best of our knowledge, been successfully nu- 
merically demonstrated in the literature previously. In 
the most extensive numerical study of the Smoluchowski 
equation to date, that of Lee j26l . the numerical difficul- 
ties posed by kernels like Eq. ([7]) were explored and it 
was concluded that there are no self-consistent solutions 
of Eq. ([T]) for such kernels. Our results demonstrate that 
the non-conservative regularisation, Eq. ([S]), provides one 
way around these difficulties, at least numerically, al- 
though we expect that the regularisation could also be 
useful for rigorous mathematical studies. 

From a physical point of view, the most important ob- 
servation about the results presented in Fig.[3]is that the 
regularised gelation time decreases extremely slowly as 
the cut-off is increased. t*]yj decreases by a factor of less 
than 2 as M is increased by 8 orders of magnitude. This 
very weak dependence means that, in practice, kernels 
which would exhibit instantaneous gelation, even com- 
plete instantateous gelation, for M = cx), can still have 
smooth, physically reasonable behaviour for finite M. 
The regularised gelation time, t^, depends sufficiently 
weakly on the actual value of M that such regularised 
systems may still be useful in modelling. On the basis 
of our numerics, we conjecture that as M — >■ oo, the reg- 
ularised gelation time decreases as ~ (logM)~" for 
some a > 0. This would complement heuristic arguments 
put forward by Ben-Naim and Krapivsky |30j for gela- 
tion in finite systems of particles undergoing exchange- 
driven growth for which it is argued that the gelation 
time decreases as a power of the logarithm of the initial 
number of particles when the aggregation rate increases 
sufficiently quickly as a function of the mass of the larger 
cluster and thus becomes instantaneous in the limit of an 
infinite number of initial particles. An important ques- 
tion, which is left open, is to develop an analytic approach 
allowing the determination of the functional dependence 
of ^A/j^^,^^ on Mmax and the value of a if the dependence 
conjectured above in indeed correct. We feel it is un- 
likely that the numerics can be extended to sufficiently 
large values of M^ax to determine this dependence un- 
ambiguously although on the basis of the numerics we 
have available and taking into account the analytic work 
on the corresponding problem with a source of monomers 
reported below, we would not be surprised if a = 1/2. 
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FIG. 4. (Color online) Total particle density, Mo{t) 
/q°° Nm{t) dm for K(m\,m2) = mi + m2 and source. 



We devote much of the remainder of the paper to 
studing what happens when a source of monomers is in- 
troduced into a system with an instantaneously gelling 
kernel. Such a situation has been partially analysed by 
Kontorovich [T7| in the context of gravitational cluster- 
ing and by Horvai et al. [H] in the context of differential 
sedimentation. Both studies concluded that the system 
should reach a stationary state at large times in which 
injection of monomers is balanced by the aggregation of 
smaller clusters into larger ones (a mass "cascade"). In 
both cases, the cascade is non-local in the sense that the 
transfer of mass from small clusters to large is dominated 
by the interaction of very large and very small clusters. 
For a detailed discussion of the criteria for locality of 
mass transfer in cascade solutions of the Smoluchowski 
equation see [3T] . 



III. THE ADDITION MODEL: A SIMPLIFIED 
DESCRIPTION OF RUNAWAY ABSORPTION 
OF SMALL CLUSTERS BY LARGE ONES 



Instantaneous gelation is driven by the runaway ab- 
sorption of small clusters by large ones. This fact is 
most easily seen from the analytically tractable (but 
non-gelling) marginal kernel, mi -f m2, with source of 
monomers. For e = 0, the kernel in Eq.([7| is the stan- 
dard sum kernel. In this case, integrating Eq.([T]) with 
respect to m (we have set toq = 1) allows us to obtain 
the following equation for the total number of particles 
in the system: 



dN 
~dt 



There is no gelation transition at finite time in this case 
so Aii{t) = Jt and we get a closed equation for N{t): 



dN 



-2JtN{t) + J. 



The solution is 



(9) 



(10) 



where Erfi(a;) = /o™ *^ dy is the "imaginary er- 

ror function" . This formula shows an interesting feature 
which is surprising at first sight: while the total num- 
ber of particles in the system initially increases as we 
add particles, it subsequently reaches a maximum and 
starts to decrease as shown in Fig. |4] It tends to zero 
(x 1/t) as time gets large. Meanwhile the total mass 
in the system increases linearly. This tells us that the 
long time behaviour of the system is dominated by an 
ever decreasing number of increasingly massive particles 
which immediately eat all the monomers injected into 
the system. Thus we already see the essential feature of 
non-local interactions at the level of the sum kernel. Fig. 
m demonstrates this behaviour. The numerical solution 



follows the analytic prediction, Eq. (10), for longer and 



longer times as the mass cut-off M is increased. This case 
is marginal in the sense that there is no finite time gela- 
tion even though all the mass gets concentrated in larger 
and larger clusters. If e > 0, and the exponent v > 1, 
then big clusters are so "hungry" that they consume all 
smaller clusters at a rate which diverges with the cut-off. 

The dynamics of a coagulation system dominated by 
the absorption of small clusters by large ones is modelled 
in extremis by the so-called "addition model" (Brilliantov 
& Krapivsky |32) ) in which clusters are only permitted 
to grow by reacting with monomers. Many of the inter- 
esting phenomena which we observe in the full coagula- 
tion problem have qualitative analogues for the addition 
model so we shall devote some time to studying it. The 
addition model corresponds to the kernel 



K{mi,'ni2) 



^ (k(toi) S{m2 



1) + K{ni2) S{mi 



1)), 
(11) 

where the function ^(m) is typically taken to be a homo- 
geneous function of degree v. This considerably simplifies 
the coagulation dynamics. Although the addition model 
can be analyzed for a fairly general set of gelling reaction 
rates {K(m)},„>i, we will use the following example for 
illustrative purposes: 



K,{m) 



1 



m 
m 



1 

1,2, 



, M, e > 0. 



(12) 



-2Mi{t)N{t) + J. 



(8) 



Here we have changed notation slightly and written v = 
1 -|- e so that the instantancnous gelation regime corre- 
sponds to e > 0. The reason for introducing a separate 
parameter, 7, to describe the monomer reaction rate will 
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become clear in Sec. |VI| where we show that by choos- 
ing 7 appropriately, it is possible to make the connection 
between the addition model and the Smoluchowski equa- 
tion more quantitative. We shall also impose the initial 
condition 



(13) 



excluding the possibility of zero initial momomer con- 
centration (this is to regularize a change of time to be 
defined below). 

Insertion of Eqs. (Ill and (12 1 into Eq. ([!]) yields the 
discrete set of equations 



dNi 
dt 



M 



it) = -2 7 _ ^ «(m)iViiV,„ + J, (14) 



(t) = K{m - l)NiNjn-i - K{m)NiN„ 
2 < TO < M, 



(15) 



where explicit dependence of the on t has been sup- 
pressed for brevity and a cut-off, M, has been introduced. 
Instantaneous gelation in the addition model without a 
cut-off for z/ > 1 in the absence of a source was conjec- 
tured in [21] and subsequently proven by Laurencot [33] 
providing support for our intuition that the process of 
aggregation of monomers by large particles is a runaway 
process in the full aggregation equation when ly > I. Be- 
fore continuing, we would like to emphasise again that 
we only expect the addition model to be similar to the 
Smoluchowski dynamics in the instantaneously gelling 
regime when the dynamics is dominated by the aggre- 
gation of large clusters and monomers. It is known that 
the addition model has very different dynamics to the 
Smoluchowski equation in the non-gelling regime, v < 1. 
In particular it does not exhibit scaling in the absence of 
a source [5^ . 

Eqs. (fT4| 15 1 simplify after a non-linear change of 



physical time to "monomer" time t defined as follows: 



dT{t) = Ni{t)dt, 
r(0) = 0. 



(16) 
(17) 



Once the solution parameterized by the monomer time 
has been found, the inverse map to the physical time is 
given by 



t{r) = 



dr' 



iVi(r') 



(18) 



This integral is proper at r = due to the initial condi- 
tion, Eq. ( jlS) ). In terms of the monomer time, Eqs. (14 
151 take the following form: 



M 



N[{T)^^2-iN,{r)-Y,K{m)Nm + ^. (19) 



N'^(t) = K{m - l)iV„_i - K{m)N„ 
2 < TO < Af , 



(20) 



where ' denotes the monomer time derivative Note 



that Eq. ( 20 ) is an inhomogeneous linear system of ODE's 



with the inhomogeneity being a nonlinear function of the 
monomer concentration Ni{t). It can be solved with re- 
spect to the concentrations of polymers using the Laplace 
transform in the monomer time. The solution can be 
written in the form 

poo 

N,n{T) = / dT'K„AT-T')Nl{T'), 



Jo 

TO : 



1,2,. 



where 



dbJ 



m ^ 

fc=2 



(21) 



(22) 



the integration contour F = (ct — ioo, a + ioo), and a is 
a real constant: a > — min,„>2[K(TO)]. It is easy to see 
from Eq. ^ that 



K{t) = 0, T < 0, 

which is an expression of causality in the addition model. 
Substituting ( 21 ) in ( 19 1 we arrive at a nonlinear integral 



equation for the monomer concentration: 

iV{(r) = -2 7iVi(T)-«(l) / dT'QMiT-T')N,{T') 

Jq 

J 



where 



M 

Qm{t) - J2 



duj 
— f 
2tt 



rJJ 1 

k=2 



(23) 



(24) 



It is clear from Eq. (23) that Ni{t) stays positive at all 



times. The transformation to monomer time is therefore 
well defined. To find the steady state solution of Eq, 
we note that all the poles in the integrand of Eq. 
have negative real parts. Thus, 

/■oo 

lim I dT'QM{T-T')Ni{T') = Ni{oo) j drQMir). 



Substituting this answer into Eq. 
iV'(cx)) = we find that 



23]) and setting 



J 



j{2 + j;^QM{r)dr) 



J 



7 (M + 1) ' 



(25) 



where the last equality results from direct integration of 
the kernel, Eq. (24 1, over time using Cauchy's theorem. 



For the specific case of the kernel, Eq. (12 1, the corre- 



sponding steady state answer for the polymer concentra- 
tions is 



A^,„(c3o) =77Vi(cx))k(to)~i 

-l-e 



J 



M + 1 



(26) 
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We acknowledge that it would be easier to find the steady 
state directly from Eqs. (ImI [T5|, but Eq. p3|) wiU be 



better suited to study the time evolution in what follows. 



IV. THE FULL COAGULATION PROBLEM 
WITH > 1 IN THE PRESENCE OF A SOURCE 
OF MONOMERS 

In the presence of a source of monomers, Eq. ([T]) has a 
formal stationary solution which scales for large masses 
as Nm ^ m ~ which describes a cascade of mass from 
small masses to large jSTj. This solution is only valid if 
the collision integral is convergent. For the case of the 
kernel Eq. Q, this convergence criterion fails for e > 
[3T] presumably reflecting the tendency for the system to 
gel instantaneously in this regime. When a cut-off is in- 
troduced, a stationary state may be reached if a source of 
monomers is present although this stationary state must 
involve the cut-off rather than being of the cascade type. 
In such a stationary state was found for the case of 
differential sedimentation (j/ = |) which scaled for large 

_ 4 

masses as N„i ^ m 3 . 

In this section we perform a systematic analysis of the 
stationary state for the model kernel 



K(mi,m2) = max(TOi, , 



(27) 



which captures the essential features of instantaneously 
gelling systems while permitting some convenient simpli- 
fications of the equations. If we assume that the system 
reaches a stationary state and define 



M 



r 



M 



(28) 



then some manipulations of the discrete analogue of 
Eq. ([T]) show that the stationary state can be expressed 
via a recursion relation: 



N„ 



\ Y."p=i '^ax(p, m - pYNpNm-p 



(29) 



which can be iterated to find Nm once one observes that 
the monomer concentration is fixed in the stationary 
state as: 



J 



(30) 



For each value of M, this iteration procedure can, in 
principle, be used to self-consistently determine Tm- In 
this paper, we adopt a different approach which is based 
on the assumption that the mass transfer is nonlocal in 
mass space. Following the approach of [T5j, Eq. ([T]) can 
be approximated by: 



dN„ 
dt 



M,N„ 



(31) 



where the moments A4i and are now to be computed 
with the cut-off retained: 



M 



m" Nm dm. 



(32) 



The source and sink terms have been omitted for now. In 



obtaining Eq. (31 1 it has been assumed that the integrals 



A^i and Mi, are dominated by their lower and upper 
limits of integration respectively. The consistency of this 
assumption will be determined a-posteriori. Since A^i 
and j\4i, do not depend on m, Eq. (31 1 has the stationary 
solution: 



C exp 



/3- 



,1-1/ 



1 



TO 



(33) 



where C is an arbitrary constant and we have, for con- 
venience, introduced the parameter. 



Ml' 



(34) 



Although we do not know the value of j3 a-priori, it can 
now be determined self-consistently due to the fact that, 
when Nm is given by Eq.(33), the integrals M.i and M.u 



can be expressed explicitly in terms of incomplete gamma 
functions: 



Mm 

MAP) 



C fl-v 
C (\-v 



1 PM^-"" 



V - 1 



-r 



V 



v-l v-l 

_l ^ 



After some algebra, Eq. (34) reduces to the following consistency condition for the value of /3: 
-2 ^Afi-"^ 



I/- 1' v-l 



v-2 /3 
v-r^v-l 



I 

v-l 



13 



I 



l-v' v-l 



I 

v-l 



l-v' v-l 



(35) 
(36) 

(37) 



Given v > I, this can be solved numerically for any value of M. Such a numerical investigation indicates that j5 is 
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a slowly increasing function of M. Furthermore, since 
we expect the integrals M.i and J^i, to be dominated 
by their respective lower and upper limits of integration 
as M — T' oo, we would expect Eq. (37 1 to satisfy the 



following asymptotic balance for large M : 



1' 



1 



1 



1 



1 



(38) 

Since numerics indicate that (3 increases (but only slowly) 
as Af grows and > 1, the argument of the left-hand 
gamma function, should increase in the limit of in- 
terest. Likewise, the argument of the right-hand gamma 
hmction, ^ , should decrease in the limit of interest. 
The relevant leading order asymptotics are then: 

r(a,z) - z"-i e-^as z -> oo (39) 
T{a,z) i z° as z -> 0. (40) 



Substituting these into Eq. (381, one finds that the lead- 



ing order terms balance as M — )■ oo provided 
13 ^\og{Ar~^). 



(41) 



Whilst the index v — 1 could be absorbed within the 
implied front factor in /?, we keep it in the following ex- 
pressions since it captures how the dependence on the 
cut-off, M, goes away when v approaches unity. Know- 
ing /3, Eqs. (40) and (391 can now be used to detemine 



the asymptotic behaviour oi Aii and M.i, from Eqs. (35) 
and (36 1. After some work one finds: 



Ml - 



CM 



\og(M'^- 

- CM 



|-as M — >■ oo 
as M oo. 



(42) 
(43) 



It now remains to determine the constant C in Eq. ( 33 ) 



in order to close the argument and allow us to check for 
consistency. This can be done using global mass balance. 
Multiplying Eq. ^ (with the kernel Eq. ([27])) by m and 
integrating from 1 to M yields the global mass balance 
condition: 



J = 



M 



dm m N„ 



M 



dmi miN„ 



(44) 



Our assumption of nonlocality allows to replace the inner 
integral by in the limit of large AI since it is domi- 
nated by its upper limit. Consequently the global mass 
balance condition in the limit of large M is 



J = MiM^ = 



log(M'^-i) 



where we have used Eqs. (42) and (43). This gives 



C 



VJlog(Af'--i) 
M 



(45) 



(46) 



Putting this all together with Eq. (33) the asymptotic 
solution of Eq. (31 ) is 



(47) 




10" 10 



10^ 10^ lO"* 10^ 10^ 10^ 
m 



FIG. 5. (Color online) Non-local stationary state (theory vs 
numerics) for v — 3/2. 



Everything becomes self-consistent: A4i becomes inde- 
pendent of AI and Aii, diverges as y/log{AI'^~^) when 
AI gets large thereby justifying the use of Eq. (jsTJ. This 
theory does a good job of capturing the general features 
of the stationary state: stretched exponential decay at 
small cluster sizes followed by a cross-over to a power law 
decay with exponent i' for large cluster sizes but with an 
amplitude which decreases with the cut-off, M. Fig. [5] 
compares Eq. (47) with the stationary state obtained 



from numerical simulations of Eq. ^ with ly = 3/2 for a 
range of values of AI. The agreement is excellent given 
that there are no adjustable parameters in Eq. (47). Note 
that the large mass behaviour, N.^ ^ m"", is in agree- 
ment with the stationary solution of the addition model. 



Eq. (26), obtained in Sec. Ill 



Finally, we can now quantify the rate at which the 
system becomes singular as the cut-off, Af , is increased. 
Using Eq. (47) to calculate the total particle density in 



the stationary state we get 



Af- A//^'^' 



N 



Vj(^A'I 

Af yiog(M^ 




as Af — oo. 

(48) 

The density of particles in the stationary state thus van- 
ishes as the cut-off is removed which is the signature of 
instantaneous gelation. The approach to zero is logarith- 
mically slow however. The asymptotic estimate, Eq. (48 1, 



is compared against the total density measured from a 
sequence of numerical simulations with v = 3/2 and dif- 
ferent values of the cut-off in the inset of Fig. [5] The the- 
ory again produces convincing agreement with numerics 
without any adjustable parameters. 
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FIG. 6. (Color online) Total density vs time for K{mi, 7712) 



O 




FIG. 7. (Color online) "Q-factor" of the oscillatory transient 
for a range of values of 1/ = 1 + e and a range of cut-offs, M. 



with nucleation and coagulation of droplets [331 US] which 
may be related to the effect noted here. From Fig.|6j it is 
not obvious that these oscillations are decaying at all for 
e > 0. From our numerical investigations we believe that 
these oscillations are indeed transient although the rate 
of decay becomes very slow as the cut-off is increased. 
Although the oscillations are quite nonlinear in charac- 
ter, we defined a quantity analogous to the Q-factor of 
a linear oscillator by measuring the ratio of subsequent 
maxima of the signal and defining the Q-factor to be the 
estimated numerical limit of this ratio as t — > 00 (for all 
of our numerical experiments with e > 0, we found that 
this ratio does indeed seem to converge to a finite value) . 
Thus a Q-factor of 1 would correspond to persistent os- 
cillations. 

A summary of these numerical experiments is pre- 
sented in Fig. [7] which shows numerically estimated Q- 
factors for the oscillations observed for a range of values 
of e and M. The most evident trend from the data is 
the fact that the oscillations come closer and closer to 
a Q-factor of 1 as M is increased. These measurements 
indicate that the oscillations are decaying in all cases but 
only very slowly when the cut-off becomes large. The be- 
haviour as a function of e is less clear with some possible 
evidence for a weak maximum. 

We offer the following heuristic explanation for how os- 
cillations can be generated in this system. For e > 0, the 
case for which instantaneous gelation should occur in the 
absence of a cut-off, the indication of Fig. [3] is that the 
particle density drops close to zero in finite time (but this 
time is not zero). When monomers are injected into the 
system, until the cut-off is felt, the mass in the system 
just increases linearly. However, as soon as large particles 
are generated, their absorption of the monomers is a run- 
away process which very rapidly converts the monomers 
which have accumulated in the system into large particles 
which are immediately removed by the cut-off. This then 
resets the system close to its initial state in which there 
are almost no particles in the system. The dynamics then 
repeats. 



V. APPROACH TO THE STATIONARY STATE 

Numerical simulations indicate that the approach to 
the stationary state is interesting. Fig. |6] shows a plot of 
the total density in the system as a function of time for 
the generalised sum kernel, Eq. ([t]), for a range of val- 
ues of e with cut-off AI = 10^. In the nonlocal regime, 
e > 0, the approach to the stationary state is charac- 
terised by long-lived transient oscillations of the density 
in the system. This is surprising since there is nothing of 
an oscillatory or cyclic character in the underlying coagu- 
lation dynamics. We note, however, that long-lived oscil- 
lations of the total droplet density have been observed in 
experiments on the phase separation of binary mixtures 



VI. APPROACH TO THE STATIONARY STATE 
IN THE ADDITION MODEL 

We do not, as yet, have an analytic characterisation 
of the oscillations described in the previous section. In- 
stead, we conclude this article with an analysis of the 
analogous problem for the simpler case of the addition 
model discussed in Sec. |III| for which the oscillatory ap- 
proach to the stationary state, Eq. (26), can be derived. 



The addition model aims to model the most impor- 
tant contribution to the dynamics of the full aggrega- 
tion model with non-local kernels by concentrating on 
interactions between "light" and "heavy" particles only. 
At the level of modeling, there is no reason to expect 
that the parameters of the addition model are the same 
as the analogous parameters in the full model. Before 
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proceeding with calculations, we first present a semi- 
heuristic argument for the most appropriate choice of 
parameters. Let us first clarify what is meant by "light" 
and "heavy" in terms of the steady state of the full 
model, Eq. (47). Heavy particles have masses in the 
polynomial tail of Eq. (47|, m > log(A/)''"^ Light 



particles belong to the stretched exponential region of 
Eq. (47), m <^ log(M)'^~^. Therefore, neglecting log- 
arithmic corrections, we can treat all light particles as 
effective monomers. Hence, the parameters J and Km 
for TO > 1 in the addition model are the same as in the 
full model, but the parameter 7 — k{1) should be cho- 
sen as an effective interaction rate of all light particles 
in the full model with themselves. The easiest way to 
choose this rate is to require that the steady state of the 
addition model gives the correct description of the mass 
distribution of the full model at large masses. Compar- 



ing Eq. ( 26 ) with the large-TO behaviour of Eq. ( 47 ) we 
see that we should choose: 



7 ' 



1 

M' 



(49) 



For large times let us look for a solution of Eq. ( 23 ) in 
the form 



iVi(r) 



J 



7 {M + 1) 



{1 + Ae<^), 



(50) 



where ReC < 0. The limit of large times formally corre- 
sponds to the limit of small amplitude A. The require- 
ment that Eq. ( 50 ) solves Eq. ( 23 ) up to terms of order 
leads to the following non-linear eigenvalue problem: 

M 

C + 7 (M + 3) + 7 ^ ^"=^ = 0. (51) 

Notice that for gelling kernels. 



lim V log ( 1 



fc=2 



c 



K{k) 



< 00 



Therefore, the sum in the left hand side of Eq. (51) di- 
verges as the first power of M: 

1 

lim — V e~^"=2'°s(i+^) ^ j:~^iog(i+^) 

M-i-oo M ^ 
m=2 

As we are interested in the limit of large cut-off mass M, 
let us solve the eigenvalue problem, Eq. (51), in the 
vicinity of M = 00. Bearing in mind that 7 ~ ^ the 
leading term of the large M-expansion of the left hand 
side of Eq. (51 ) gives 



g-log(H-C)-Er=2 los(l+;jiE)) 

or equivalently, 

log I i -t- - ^ 

fc=2 



-1, 



(52) 



iog(i+c)+f:iog(i+^) 



= 7ri(l + 2(7) (53) 



where 5 €! Z. We need to find the solution to the above 
equation with the largest real part. It turns out that an 
analytical solution is possible for the set of reaction rates 
Eq. (12) in the limit e ^ 1. Let us Taylor expand the 



logarithms in Eq. ( 53 ) and switch the orders of summa- 



tion in anticipation of seeking an asymptotic expansion 
for C. We obtain: 

CO 

^(-l)"+i5„C" = vrz(l + 2g) 



where 



n 



fc=2 



For the reaction rates given by Eq. (12) with e <C 1, we 



have the following behaviour for the Sn 



1 



1, n = 2,3, 



The fact that 5*1 ^ 5„ for small values of e allows one 
ad the s 
or (7 



to find the solution to Eq. ( 53 ) with the largest real part 
(q — ov q = —1) as an asymptotic expansion in e: 

:S2 



C = ±t 



1 



^3 r, 3 S2 



+ (-7r^ 



(54) 
Therefore 



We observe that ReC < and |ImC| » |ReC| 
the monomer concentration slowly decays to the station- 
ary value in an oscillatory manner with many periods 
of oscillations per inverse decay rate. Accordingly the 
oscillations' quality factor is close to 1 for e ^ 1: 



(55) 

The independence of the quality factor of the cut-off mass 
M (up to possible logarithmic corrections) is confirmed. 
It is also straightforward to check that the physical period 
of oscillations is A/- independent: notice that Re^ < 0. 



Therefore, at large times the transformation Eq. ( 16 ) 
reads: 

t(0 =To + A^i(oo)t + 0(e^°^*). 



Substituting this expression into Eq. ( 50 ) we find that in 
the limit M 



■ 00, 

Ni{t) = \/j(l + Ae«^*), 



(56) 



where C is given by Eq. (54). Therefore, the period of 



oscillations does not depend on the cut-off. It scales as 
the inverse square root of the flux of monomers: 

T= -^5i(l + 0(e3)). (57) 

The period diverges in the limit e — >■ cxi, which implies 
that oscillations disappear as we increase the locality of 
the aggregation kernel. These conclusions are in agree- 
ment with our numerical measurements on the full coag- 
ulation equation presented in Sec. [Vj 
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VII. CONCLUSIONS 

To conclude, the fact that certain aspects of the regu- 
larized Smoluchowski equation, such as the total density, 
depend so weakly on the value of the cut-off, M, used to 
regularise the system means that instantaneously gelling 
kernels are still potentially reasonable as physical models. 
Their use to describe gravitational clustering or differen- 
tial sedimentation is not in contradiction with mathe- 
matical results which indicate that the density vanishes 
for all positive time in the unregularised system provided 
one is willing to accept non-universal dependencies on a 
cut-off as part of the model. We have presented an essen- 
tially complete analysis of the stationary state of the reg- 
ularised system in the presence of a source of monomers 
when the regularisation is done be removing all clusters 
having size larger than the cut-off. In this case, a power 
law scaling is observed for large cluster sizes whose expo- 
nent is universal (depending only on the scaling proper- 



ties of the kernel) but whose prefactor is strongly cut-off 
dependent. The scale-to-scale mass balance is of a non- 
local character rather than being described by a mass 
"cascade" as is the case for regular gelling and non-gelling 
kernels. Finally, we found that this stationary state is ap- 
proached in an interesting way with long-lived oscillatory 
transient resulting from interaction between the cut-off 
and the source. A quantitative analysis of the approach 
to the steady state was provided in the simplified case 
of the addition model which expresses the most extreme 
case of nonlocal mass transfer in which clusters can only 
grow by interaction with monomers. 
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